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Abstract. We investigate the relation between the dynamical features of a 
supercooled liquid and those of its potential energy landscape, focusing on a model 
liquid with density anomalies. We consider, at fixed temperature, pairs of state points 
with different density but the same diffusion constant, and find that surprisingly they 
have identical dynamical features at all length and time scales. This is shown by the 
collapse of their mean square displacements and of their self intermediate scattering 
functions at different wavevectors. We then investigate how the features of the energy 
landscape change with density, and establish that state points with equal diffusion 
constant have different landscapes. In particular, we find a correlation between the 
fraction of instantaneous normal modes connecting different energy minima and the 
diffusion constant, but unlike in other systems these two quantities are not in one-to- 
one correspondence with each other, showing that additional landscape features must 
be relevant in determining the diffusion constant. 
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Particle motion in glassy systems is frequently conceptualized as consisting of caging 
periods, during which the system rattles around in a energy minimum of its energy 
landscape, interrupted by transitions between minima. In this scenario the features 
of the potential energy landscape (PEL) play a fundamental role in determining the 
overall dynamical properties. One way to investigate how a system explores its 
PEL is via the study of the instantaneous normal modes (INM) [TJ, which are the 
instantaneous eigenvalues and eigenvectors of the Hessian matrix. The eigenvalues give 
information about the curvatures of the PEL hypersurface around the point the system 
is visiting, while the eigenvectors give the corresponding directions for the joint motion 
of particles through phase space. As concerns the dynamical features of the system, 
the negative eigenvalues are of particular interest as they identify the unstable modes 
of the system M- Specifically, some of these negative eigenvalues are expected to be 
directly related to diffusion. It is clear, however, that not all of the unstable modes 
can contribute to diffusion as unstable modes are also observed at temperatures which 
are sufficiently low for the diffusion constant to be negligible. Strategies must therefore 
be devised to determine which are the unstable modes that are relevant to particle 
diffusion [U El E] • 

Here we explore a proposal by La Nave and coworkers [5j, who suggested that the 
unstable modes allowing for diffusion are those connecting different energy minima, 
known as double-well modes. Indeed, a one-to-one relation between the diffusion 
constant and the fraction / DW of double-well modes has been observed in a range of 
model systems [3EUH], We note that a different approach, where the modes contributing 
to diffusion are identified with the unstable and extended ones, and detected via a finite- 
size scaling analysis |2j, has recently been shown to yield analogous results [6]. 

In order to investigate the relation between the features of the potential energy 
landscape and the diffusion constant of liquids, we focus on a simple model characterized 
by water-like density anomalies, i.e. regions of the phase diagram where the dynamics 
speeds up upon isothermal compression. We first consider two state points with the 
same diffusion constant D and same temperature, but different densities, and show 
that key dynamical features coincide at all length and time scales. Indeed, not only 
are the diffusion constants equal, but also the mean square displacements and the self¬ 
intermediate scattering (i.e. incoherent correlation) functions at different wavevectors 
coincide at all times. Differences are observed only in four-point quantities such as the 
non-linear dynamical susceptibility. 

We then consider how the PEL changes with density, and how these changes are 
related to those in the dynamics. We show that the typical distance between two energy 
minima connected by a double-well mode, as well as the associated typical energy 
barrier, increase monotonically with density. Both quantities must thus be different 
for two state points at different density, even if they have the same diffusion constant. 
On the other hand, the fraction of unstable normal modes connecting different energy 
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minima /dw varies non-monotonically with the density. This ensures it is correlated 
with the diffusion constant, but our results do not support the presence of a one-to-one 
relation between D and /dw such as has previously been seen in other systems. 

2. Model system 

As a simple model system for investigating the relation between the dynamics and the 
features of the energy landscape we consider a polydisperse mixture of N harmonic 
disks of identical mass m, in two dimensions. Diameters are uniformly distributed in 
the range [d min , d max ], with the difference d max — d rnin between the largest and smallest 
diameter being 82% of the mean diameter (d max + d min ) /2 so that the distribution is 
fairly broad; this is necessary to prevent crystallization. Two particles i and j interact 
via the repulsive harmonic potential 



( 1 ) 


The interaction is of finite range: particles interact only if their distance r is smaller 
than their average diameter d = (</ + </,)/2. Lengths, masses and energies are expressed 
below in units of d max , m and of e, respectively, and the density is specified in terms of 
the volume (or more precisely, area) fraction </> = N(A)/L 2 . Here L is the system size, 
(A) the average particle area, and N the number of particles. 

We have performed molecular dynamics simulations [HI] of this system with 
N = 10 3 particles, integrating the equations of motion using the Verlet algorithm, with 
a timestep of 5t = 10“ 4 . Dynamical quantities are averaged over 10 3 independent runs. 
The system is first brought to thermal equilibrium via simulations in the NVT ensemble, 
while production runs are carried out in the NVE ensemble. Thermal equilibration 
is ascertained by checking for the absence of aging. The INMs are determined by 
equilibrating a larger sample, consisting of IV = 4 x 10 3 particles, and averaging over 50 
independent realizations to obtain accurate statistics. 

3. Density anomaly 

In the majority of molecular and colloidal fluids the dynamics slows down upon 
isothermal compression. There are, however, fluids that display the opposite behavior 
in some region of their temperature-density phase diagram. Such density anomalies 
are found in a number of fluids; the most important of these is certainly water, with 
other examples including Si, Ge, Sn, and ionic melts with suitable radius ratios, such 
as Si0 2 , BeF 2 and Ge0 2 [11] [12]. While all of these systems have an anisotropic 
interaction potential that favors short-range tetrahedral order around the particle 
centers, density anomalies have been also observed in models with spherically symmetric 
interaction potentials [13]. These systems have been studied extensively to identify the 
physical origin of density anomalies and to connect them to specific features of the 
interaction potential. The harmonic potential considered in this work is an example 
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Figure 1. Non-monotonic density dependence of the diffusion constant. The arrows 
indicate the states Hi = (</> = 1.0, T = 0.18) and H 2 = (</> = 1.3, T = 0.18), which 
have similar diffusion constants. 

of a radially symmetric potential giving rise to a density anomaly, as illustrated by 
the non-monotonic dependence of the self-diffusion constant on volume fraction shown 
in Fig. [I] We have previously investigated the physical mechanisms responsible for 
this behavior, and concluded that it must be attributed to the emergence of contacts 
between a particle and those of its second coordination shell [HUS]. In the following 
we first compare the dynamical features of two state points, Hi = (4> = 1.0, T = 0.18) 
and H 2 = (4> = 1.3, T = 0.18), that have approximately the same diffusion constant, 
D H2 /D Hl ~ 1.13 (Fig. [[]), and then study the density dependence of some important 
features of the energy landscape. 

4. Dynamics 

As illustrated in Fig. [T] the presence of a density anomaly means that at a given 
temperature there are pairs of points that have different density but the same diffusion 
constant, such as state points H 1 and H 2 . Due to the difference in density, one might 
expect the dynamics at these state points to be different in general, apart from sharing 
an overall timescale set by the diffusion constant. Surprisingly, we find that the opposite 
is true: state points with the same diffusion constant share many dynamical features. 
Consider first the mean square displacement shown in Fig. [21 It is equal between the 
two state points not only at long time, where this is expected because the two states 
have equal diffusion constant, but in fact at all times to within our numerical accuracy. 
This is consistent with the existence of a close link between the Debye-Waller factor, 
which is the value of the mean square displacement in the plateau region, and the 
diffusion constant. Such a link has been put forward on more general grounds [K], and 
we have previously verified that it applies across a large part of the phase diagram 
of the model system considered here d- The mean square displacement is given 
by the variance of the self part of the van Hove function, G(r,t ), defined here as 
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Figure 2. State points Hi and H 2 have the same mean square displacement as a 
function of time (left panel), as well as the same van Hove correlation function for 
every time interval t (right panel), within numerical uncertainty. 


the probability that a particle has moved a distance r along one coordinate (say x- 
)axis at time t. To investigate the relation between the dynamics at state points H\ 
and H 2 more closely, we therefore next look at the van Hove function itself. Fig. [2 
shows that, surprisingly, the van Hove functions — evaluated for the same time interval 
t - at H\ and H 2 coincide not just in their variance, but in their entire functional 
dependence on r. The observed equality of the van Hove functions suggests that in fact 
the dynamics of the two systems might be equal at all relevant length scales. We have 
checked this by measuring the self-intermediate scattering (or incoherent correlation) 
function, F(k,t ) = j? ^X!yLi ex P [~lk • ( Vj(t ) — r*j(0))]^ for 11 different wavevectors k 
with lengths k = |fc| evenly spaced between k max = An/d max and k min = 0.27r/d max . As 
shown in Fig. [3j H 1 and H 2 have very similar self-intermediate scattering functions, 
at all times and length scales. Of course F(k,t ) is related to G(r,t), essentially by 
Fourier transform, so this equality primarily provides a consistency check that equality 
of the dynamics is found both in real and Fourier space. We note that a similar collapse 
cannot be observed in the distinct part of the van Hove function, G rf (r, t), as at t = 0 
this coincides with the radial distribution function, which is distinct for states having 
different densities. Similarly, no collapse can be observed in the total-intermediate 
scattering function. 

Digging deeper, differences in the dynamics at state points Hi and H 2 finally 
turn up in the fluctuations of the self-intermediate scattering function. If we denote 
this fluctuating version of F(k,t ) by F(k,t) = A YljLi ex P [—ik ■ ( rflt) — Tj-(O))] then 
its variance defines the so-called dynamical susceptibility x{k,t) = Ff[(F(kfl)) 2 — 
(F(k,t) 2 )]. In this quantity we do find clear differences between H\ and H 2 , as shown 
in Fig. [3l I 11 particular, while the susceptibilities attain their maximum roughly at the 
same time t* for any given k, the actual maximal value \(k,t*) is smaller for state H 2l 
i.e. in the anomalous region. 
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Figure 3. The left panel shows the self intermediate scattering function, F(k,t ), 
for k = A: m ; n + n(/c max — /c m in)/10, and n = 0,1, ..,10, from left to right. Here 
kmin = 0.27r/d max and fc max = 47 t/ d max . Data for state point Hi (full symbols) and H 2 
(open symbols) coincide at all times and wavevectors, to a very good approximation. 
The right panel shows the corresponding dynamical suspeptibilities. For clarity, we 
only show data for n = 0,5, 7, 8,9, from left ro rigth. 


5. Instantaneous normal—modes 

Above we have found strong similarities between the dynamical properties of the equal- 
diffusion constant state points H 1 and H 2 . If there is a significant correlation between 
the dynamics and the properties of the energy landscape, then this would suggest that 
the energy landscapes of states H 1 and H 2 should also have closely related features. 
More specifically this should apply to the statistics of those landscape properties that 
are directly related to the dynamics. Accordingly, by comparing different features of 
the energy landscape of the two states, we can infer to what extent they correlate with 
the dynamics. 

To characterize the energy landscape statistics we employ the instantaneous normal 
mode approach, and thus focus on the properties of the eigenvalues A; and the 
eigenvectors Su { of the Hessian matrix of the system. We consider first the distribution 
P(v) of the eigenfrequencies zy = y/\. These are shown in Fig. [H where following 
the convention in the literature we plot imaginary eigenfrequences arising from negative 
eigenvalues A* as — | zy | = — y/|Aj|. They can then be plotted on the same axis as the 
real eigenfrequencies, but remain distinct from them in such a visualization. Fig. [4] 
demonstrates clearly that the resulting two eigenvalue distributions differ between the 
state points H 1 and H 2 . This is presumably driven by the difference in density between 
the two states. 

The diffusive dynamics of the system, however, is not influenced by all modes. La 
Nave and coworkers [5] suggested that the unstable modes allowing for diffusion are 
the double-well modes, i.e. those connecting different energy minima, and observed a 
one-to-one relation between the diffusion constant D and the fraction of modes, /dw, 
that have this double-well character. In order to determine if an unstable mode i is a 
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Figure 4. Probability distribution of the eigenfrequencies of state points Hi and H 2 , 
for a system of N = 4,000 particles. Imaginary eigenfrequencies v are plotted as —\v\ 
to separate them from the real eigenfrequencies. 



Figure 5. For each unstable eigenvalue i we have displaced the particle positions by 
an amounty s 5ui along the eigenvector Su f, the eigenvector is normalized so that s 2 
gives the total squared displacement of all particles. The dependence of the energy of 
the system E on s allows one to discriminate between double-well modes and shoulder 
modes, as illustrated in the figure for two typical cases. The figure refers to state Hi 
and N = 4, 000. 


double-well mode, one asks how the interaction energy of the system E t (s) changes as 
particles are displaced by s 8u t along the eigenvector. Note that, since we consider unit 
eigenvectors, s measures the distance the system is displaced in phase space: s 2 is the 
sum of all squared particle displacements. If E^s) is a double-well function, i.e. has two 
local minima, then we classify the mode i as a double-well mode. As an illustration, 
we show in Fig. 0 numerical results for E^s ), for the two cases of a double-well and of 
a “shoulder” mode. We note that all modes with approximately zero eigenvalue, which 
are responsible for the observed peak in P{y) around v = 0, turn out to be shoulder 
modes. 

The features of the double-well modes change with density. As displayed in 
Fig. H the average distance (dr) separating two energy minima in phase space and 
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Figure 6. Bottom row: Probability distribution of the distance dr between two energy 
minima connected by a double-well mode (bottom left), and of the energy barrier A E 
separating the two minima (bottom right), for the two states H\ and H 2 . Top row: 
averages of these quantities as a function of volume fraction, across the range between 
Hi (</> = 1) and H 2 (4> = 1.3). 

the corresponding average energy barrier (A E), increase monotonically with density. 
In particular, these quantities are therefore not identical for the two states H 1 and 
H 2 , suggesting that they are at most weakly correlated with the dynamics. The figure 
(bottom row) also show that dr is approximately Gaussian distributed, while A E has a 
roughly exponential distribution. We already know the averages of the two distributions 
are different for H\ and H 2 , and accordingly also the distributions themselves are quite 
distinct. 

So far we have found for the instantaneous normal mode eigenfrequency 
distribution, and geometrical properties of double-well modes such as distance and 
barrier between the two minima, that states H 1 and H- 2 show clear differences in their 
energy landscape statistics. We now turn, finally, to the fraction f DW of double-well 
modes. As shown in Fig. [TJi, this varies non-monotonically with the volume fraction, 
much like the diffusion constant. As expected, a parametric plot of D versus few 
in Fig. [7b therefore exhibits a strong correlation between these two quantities. The 
fraction of double-well modes is therefore clearly one feature of the energy landscape 
that is closely linked to the dynamics, and indeed (Fig [7b) its values in states Hi and 
H 2 are within 23% of each other. However, Fig [7b also demonstrates that for our 
system there is no one-to-one correspondence between D and fowi in contrast to what 
has been observed in other systems 

6. Conclusions 

We have investigated a polydisperse mixture of disks interacting via a repulsive finite- 
range harmonic potential, in the volume fraction and temperature range where the 
system displays a density anomaly. We have focussed primarily on the possible existence 
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Figure 7. (a) The fraction /dw of unstable modes connecting different energy minima, 
i.e. double-well modes, varies non-monotonically with the volume fraction, (b) It is 
therefore strongly correlated with the diffusion constant, but does not determine it 
uniquely. 


of correlations between the dynamics and features of the potential energy landscape. 
We have observed striking similarities between the dynamical features of state points 
characterized by the same diffusion constant: as far as can be ascertained from 
incoherent (i.e. single particle) two-point correlations, their dynamics appear to coincide 
at all length and time scales. We note that these dynamical similarities are not confined 
to the two specific state points Hi = (</ = 1,T = 0.18) and H 2 = (4> = 1.3, T = 0.18) 
that we have studied; we have obtained analogous results when comparing e.g. the 
dynamics at the state points (</> = 1.03, T = 0.16) and (</> = 1.25, T = 0.16). What 
is the physical origin of these dynamical similarities is not obvious to us, and will be 
an interesting question to follow up in future work. In this respect, we mention also 
that not all liquids with density anomalies have this feature; indeed, for a repulsive 
Gaussian interaction we generally find quite distinct dynamics for states of equal 
diffusion constant [15]. 

Despite the dynamical similarities we have found in our specific model system, 
the states that are paired up by having equal diffusion constant have PELs with quite 
different statistics. These include the average distance between minima connected by 
unstable normal modes, and the average barrier. One quantity where we do find a strong 
correlation with the dynamics is in the fraction /dw of double-well modes among all 
instantaneous normal modes; such double-well modes are defined as unstable modes 
connecting different energy minima. However, while there is a clear correlation between 
/dw and D, the two quantities are not uniquely determined by each other. This suggests 
that in our system there must be other features of the landscape, beyond /dw, that act 
to determine the diffusion constant D. 
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